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Using an effective potential method we examine binary black holes where the individual holes 
carry spin. We trace out sequences of quasi-circular orbits and locate the innermost stable circular 
orbit as a function of spin. At large separations, the sequences of quasi-circular orbits match well with 
post-Newtonian expansions, although a clear signature of the simplifying assumption of conformal 
flatness is seen. The position of the ISCO is found to be strongly dependent on the magnitude 
of the spin on each black hole. At close separations of the holes, the effective potential method 
breaks down. In all cases where an ISCO could be determined, we found that an apparent horizon 
encompassing both holes forms for separations well inside the ISCO. Nevertheless, we argue that 
the formation of a common horizon is still associated with the breakdown of the effective potential 
method. 
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I. INTRODUCTION 

The inspiral and coalescence of binary black hole sys- 
tems is a prime target for upcoming gravitational wave 
detectors such as LIGO. Such systems will be circularized 
by the emission of gravitational waves, and will evolve 
through a quasi-equilibrium sequence of circular orbits. 
At the innermost stable circular orbit (ISCO) we expect 
a transition to a dynamically plunging orbit. It is an- 
ticipated that this transition will impart a characteristic 
signature on the gravitational waveform. It is therefore 
important to know the orbital frequency at the ISCO, 
since the corresponding gravitational wave frequency is 
predominantly just twice this frequency. 

Predicting the waveform in detail from the transition 
at the ISCO to the final merger requires the full machin- 
ery of numerical relativity. These calculations require 
appropriate initial data. Out of the large space of solu- 
tions of the initial-value equations of general relativity, 
we need an algorithm to select solutions corresponding 
to black holes in quasi-circular orbits. The effective po- 
tential method J|] allows one to construct such solutions, 
and to determine the properties of the ISCO. 

The effective potential is based on the fact that mini- 
mizing the energy of a system yields an equilibrium solu- 
tion. This follows from the Hamiltonian equations of mo- 
tion: If the Hamiltonian TL is minimized with respect to 
a coordinate q and a momentum p, then q = dTi /dp = 
and p = —dH jdq = 0. The energy of two objects in orbit 
about each other can be lowered by placing the objects at 
rest at their center of mass. Therefore minimizing the en- 
ergy with respect to all coordinates and momenta will not 
yield a circular orbit. To find circular orbits in Newtonian 
gravity, one can minimize the energy while holding the 
angular momentum constant. This procedure works as 
well for a test-mass orbiting a Schwarzschild black hole, 



where one minimizes the ADM energy. This can be seen 
as follows. For geodesic motion, one finds 0] 




(1) 



Here M is the mass of the black hole, E is the energy 
per unit rest mass of the test-particle as seen from in- 
finity and L its orbital angular momentum per unit rest 
mass. Denote the rest-mass of the test-particle by M'. 
Then the ADM energy is simply E ADM = M + EM', and 
minimizing E ADM is equivalent to minimizing E. Hence 
minimizing the left hand side of (|l|) with respect to r 
yields the radius of circular orbits as a function of an- 
gular momentum. Minimization of (0) with respect to r 
yields r = 0, which is necessary for a circular orbit. From 
the minimum one finds the energy of the test-particle as a 
function of angular momentum. Obviously, one needs to 
keep M and M' constant during the minimization, so the 
prescription to compute circular orbits becomes: Min- 
imize -E A dm while keeping the angular momentum and 
the rest masses constant. 

These ideas have been formalized as variational prin- 
ciples for finding equilibria for rotating and binary stars 
in Newtonian gravity. There is also a similar variational 
principle for rotating stars in general relativity |». Bi- 
nary systems in general relativity are not strictly in equi- 
librium because they emit gravitational waves. However, 
for orbits outside the innermost stable circular orbit, the 
gravitational radiation reaction time scale is much longer 
than the orbital period. It is therefore a good approxi- 
mation to treat the binary as an equilibrium system. 

In this paper we apply this minimization principle to 
rotating binary black hole systems. Let the masses of 
the holes be Mi and M2, the spins be Si and S2, and 
the total angular momentum of the system be J. We 
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exploit the invariance under rescaling of the mass by 
using dimensionless quantities Mi/M 2 , Si/Mf, S2/MI, 
and J/fim, where m = Mi + M 2 denotes the total 
mass and fi — MiMi/m the reduced mass. Then we 
adopt the following straightforward prescription to locate 
quasi-circular orbits: Minimize the scaled ADM energy 
£-adm/w with respect to the separation of the holes, while 
keeping Mi/M 2 , Si/Mf, S 2 /M|, and 3/fxm constant. 

It is somewhat involved to carry out this simple pre- 
scription. The computation of the ADM energy becomes 
more difficult than for the Schwarzschild example above. 
More importantly, however, no rigorous definitions ex- 
ist for the mass or spin of an individual black hole in a 
spacetime containing two black holes. We will address 
these issues in Sec. O. Ultimately, we must use numeri- 
cal methods to generate and search among the solutions. 
Our numerical approach involves rootfinding, which is 
also described in Sec. [n| 



In Sec. Ill we present the results of the effective po- 
tential method. For the interpretation of these results, 
we need to search for common apparent horizons in our 
binary black hole data sets. These results are included 
in Sec. III. too. We discuss our results and conclusions 



in Sees. IV and [y| The appendix contains details of the 
apparent horizon searches. 



II. IMPLEMENTATION 



In order to minimize the ADM energy while keeping 
Mi/M 2 , J/fim, Si/Mj 2 and S 2 /Mf constant, we need a 
method to compute the ADM energy as a function of an- 
gular momentum, masses and spins of the holes and sepa- 
ration. As a first step we construct initial data (jij, K^) 
on a hypersurface as described in . Our particular 

approach assumes conformal flatness of the 3- metric 7^ , 
maximal embedding of the hypersurface, as well as inver- 
sion symmetry conditions on the 3-metric 7^ and on the 
extrinsic curvature Kij. The effective potential method 
is independent of these assumptions and works with all 
methods that compute initial data. For example, in Q, 
the effective potential method was used without assum- 
ing inversion symmetry. In particular, the assumptions of 
maximal embedding and conformal flatness are not essen- 
tial but merely convenient — maximal embedding decou- 
ples the Hamiltonian and momentum constraints within 
the initial-data formalism we use, and conformal flatness 
allows for an analytic solution of the momentum con- 
straints. One disadvantage of conformal flatness is that 
Kerr black holes do not admit conformally flat 3-metrics, 
at least for the simple time slicings we are aware of. In 
jjj it was shown that the Kerr metric is not conformally 
flat at second order in the spin parameter S/M 2 . Indeed, 
in Sec. Ill A we identify this deviation in our results. 

Because we assume that the initial hypersurface is 
maximal, the momentum and Hamiltonian constraints 
decouple. We follow the Bowen and York M prescription 



to solve the momentum constraint analytically. Then we 
need only solve one three-dimensional quasi-linear ellip- 
tic differential equation, the Hamiltonian constraint. It 
is solved on a so-called Cadez grid using a multigrid al- 
gorithm H . The constructed data sets depend on several 
input parameters, namely the radii and the positions of 
the throats of the holes in the flat background space, <Zj 
and Cj, i = 1,2, respectively, and their linear momenta 
and spins, Pi and S,*, i = 1,2, respectively. We note 
that in this initial-data prescription, Pj and Si represent 
the physical linear and angular momentum of the black 
hole if it is isolated. We work in the zero momentum 
frame, where P 2 = —Pi, and choose P^ perpendicular 
to C 2 — Ci in order to realize a circular orbit. Then 
the magnitude P = P\ = P 2 is sufficient to describe the 
linear momenta. Choosing a\ as the fundamental length 
scale, we are left with the following dimensionless input 
parameters: the ratio of the throat radii a — di/a 2 , the 
dimensionless background separation f3 — \Ci — C 2 |/ai, 
and the dimensionless linear momentum and spins, P/eti 
and Si/a 2 , i = 1,2, respectively. 

From the initial data we can rigorously compute the 
ADM energy £ ADM , the total angular momentum J and 
the proper separation between the apparent horizons of 
each hole, I. The total angular momentum is evaluated 
as in Ref. 0: 

J = (C 1 -0)xPi + (C 2 -0)xP 2 + Si + S 2 . (2) 

Here O represents the point about which the angular 
momentum is defined; it drops out immediately because 
Pi = — P 2 . When orbiting black holes have spin, neither 
the individual spins of the holes nor their orbital angular 
momentum L are rigorously defined. We simply take L 
to be defined by 



L = J - Si - S 2 



(3) 



with Si and S 2 defining the individual spins. 

Finally, we need to define the masses of the individual 
holes. As in Ref. [jjj, we define the mass of each hole via 
the Christoudoulou formula: 



Mf = Mf, t 
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(4) 



(5) 



where Ai is the area of the event horizon of the hole. 
Clearly this definition is only rigorous for a stationary 
spacetime. Moreover, we cannot locate the event horizon 
from the initial data slice alone. Therefore we must re- 
sort to using the apparent horizons areas in equations (^) 
and (^) instead. Apparent horizons can be determined 
from initial data and in the present case their positions 
are known to coincide with the throats of the holes [Q . 
For a stationary spacetime, apparent horizons and event 
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horizons coincide, and in a general, well-behaved space- 
time, the event horizon must coincide with or lie outside 
of the apparent horizon. In the latter case we will under- 
estimate the mass of the black hole by using the apparent 
horizon area. Some of the results of this work indicate 
that this happens for very small separations of the holes. 

With the individual masses we can finally define the 
effective potential as the non-dimensional binding energy 
of the system: 



Eb 



(E A 



Mi - M 2 )/fi. 



(6) 



Since the mass-ratio M±/M 2 is kept constant during the 
minimization, minimizing Ef,/fi is equivalent to minimiz- 
ing E AUM /m. 

We construct initial data sets starting from the in- 
put parameters a, /?, P/ai and Sj/af, and compute the 
physical parameters Eb//i, M1/M2, J / lira and Si/Mf. In 
order to construct an initial data set with certain phys- 
ical parameters we have to choose the input parameters 
appropriately. This requires nonlinear rootfinding. 

Within our effective potential approach, we will search 
for minima in the binding energy as a function of the sep- 
aration of the black holes. Fortunately, it is not necessary 
to solve for a specific proper separation t/m. It is suffi- 
cient to keep (3 constant during rootfinding and thus find 
a binary black hole configuration with some separation 
t/m. Our goal is to solve the following set of equations 
[cf. Eqns. (lOa-d) of Ref. H]: 
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The bracketed quantities on the right hand sides of (f7aj- 
[7c]) denote the physical values to be reached, and the 
expressions on the left-hand side represent functions of 
the background parameters a, P/ai, Si/af and S 2 /a\ as 
well as the fixed (3. 

For non-rotating holes, equations (7b) and ( |7c| ) are 
trivially satisfied by Si = S 2 = 0. For spinning holes 
this is no longer the case. Hence, it seems one has to 
solve the complete set of Eqns. (7a -73). However, in any 
initial data scheme where the physical spins of the black 
holes are directly parameterized, Eqns. ( |7b[ ) and ( f7c| ) can 
be eliminated. First, we note again that if the physical 
spins are directly parameterized, from Eqn. (||) we find 
that we can replace rootfinding in J /fim by rootfinding 
in L/fim. Thus Eqn. (7d) is replaced by 
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(8) 



In the zero momentum frame, Eqns. 
to 

af <xi 
Thus we can rewrite Si as 



Si 



Si Mi MiM 2 
Ml' ~M 2 L 



0- 
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and (||) simplify 
(9) 

(10) 



For a solution of Eqns. (|7a| -|7c|Jq) , the first three terms on 
the right hand side of ( |10| ) take the values of the desired 
physical parameters, so we can replace them by these 
parameters throughout the rootfinding. A similar result 
holds for S 2 ■ We perform only two-dimensional rootfind- 
ing, in a and P/a±, and set in each iteration 
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For an important subset of spin configurations, even 
one-dimensional rootfinding is sufficient as can be seen 
as follows: Consider equal-sized holes with equal spin 
magnitudes on both holes. If both spins are parallel to 
the orbital angular momentum, or both spins are antipar- 
allel, there exists a symmetry under exchange of the two 
holes. Therefore a must be equal to 1 and we are left 
with one free parameter, P/ai. If one spin is parallel 
to the orbital angular momentum and the other spin is 
antiparallel, however, this property is lost. One hole is 
co-rotating with the orbital motion and the other hole 
is counter-rotating. The choice a = 1 would result in 
holes with slightly different masses. We thus need two- 
dimensional rootfinding in a and P/ai for this case. 

Each "function evaluation" for the rootfinding involves 
the computation of an initial data set (7^ , K \j ) . High 
resolution solutions take between 30 minutes and sev- 
eral hours of CPU time on one RS6000 processor. For 
maximum efficiency, we first perform rootfinding with 
a Newton-Raphson method on low resolution data 
sets. The numerical values for Mi/M 2 and J/ fim dif- 
fer slightly between low resolution and high resolution 
solutions, therefore we solve on low resolution for ad- 
justed values of [Mi/M 2 ] and [J/ fim]. With the input 
parameters found in the low resolution rootfinding, a high 
resolution computation is performed to verify that equa- 
tions (|7^) and (7d) are indeed satisfied at high resolution, 
and to adjust the offset used in the next low resolution 
rootfinding. If necessary, this procedure is repeated. On 
average each complete rootfinding takes fewer than two 
high resolution computations. 

Following our prescription, we now minimize the bind- 
ing energy with respect to separation while keeping 
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Mi/M 2 , L/irni and Sj/M 2 constant. The binding en- 
ergy of a sequence of solutions with these quantities held 
constant represents a contour of the effective potential. 
Our code starts at large separation (3 and reduces (3 until 
a minimum in Eb/ is bracketed. Then the minimum is 
located with Brent's method H, yielding a quasi-circular 
orbit for the prescribed values of J/fim, M1/M2, and 
Si/Mf. Note that each computation of Eb/ '/j, during the 
minimization along an effective potential contour requires 
rootfinding. 

By computing quasi-circular orbits for different J / fim, 
but fixed M1/M2 and Si/Mf, a sequence of quasi-circular 
orbits is obtained. A binary black hole that radiates 
away energy and angular momentum will follow such a 
sequence approximately, assuming that the spin on each 
hole remains constant. We step towards smaller J/fim, 
and compute only as many points along each effective 
potential contour as are required for the minimization. 
As soon as we do not find a minimum in the effective 
potential contours anymore we expect to be beyond the 
innermost stable circular orbit. We trace out some com- 
plete effective potential contours around the last value of 
J / fim to check the behavior of these curves. 

Finally from the binding energy Eb/ /J. and the angular 
momentum J/fxm along the sequence, we compute the 
orbital angular frequency as 



dEb 
dJ 



The "-I — sequences" with one co-rotating and one 
counter-rotating hole. 



• The " 

holes. 



sequences" with two counter-rotating 



(12) 



We restrict ourselves to equal mass holes, M\ = M 2 = 
M with equal spin magnitudes Si = S2 = S. As we 
will see, the assumption of conformal flatness becomes 
questionable at high spins, so we consider only spin mag- 
nitudes S/M 2 < 0.50. We denote a spin configuration 
by two plus or minus signs together with a number spec- 
ifying the spin magnitude on the holes. Thus "++0.25" 
denotes a configuration with two co-rotating holes and 
spin magnitudes S1/M 2 = S 2 /M 2 = 0.25. 

Quasi-circular orbits were computed for various val- 
ues of J I fim along each sequence. In Fig. |] the binding 
energy Eb/fi along each sequence is plotted as a func- 
tion of the angular momentum J / fim. A binary black 
hole that loses energy and angular momentum through 
gravitational radiation moves along such a sequence if the 
spins of the individual holes remain constant. The dashed 
lines in Fig. [j] represent the resul ts of ( post) 2 -Newtonian 
theory which we describe in Sec. HI A . 

Using equation ([Et) we compute the orbital angular 
frequency. In Figs, g and the binding energy and the 
angular momentum along the sequences are plotted as a 
function of orbital frequency. 



sequence 



III. RESULTS 



A. Behavior at large separations 



The parameter space of spinning binary black holes is 
large - one can vary the mass ratio of the holes as well 
as spin directions and magnitudes. Astrophysically most 
interesting are holes that co-rotate with the orbital mo- 
tion, i.e. with both spins Sj parallel to the orbital angular 
momentum L. In addition to these co-rotating configu- 
rations, we examine configurations with one co-rotating 
hole and one counter-rotating hole, and configurations 
with two counter-rotating holes. We have the following 
three families of sequences: 



We compare our results to the (post) -Newtonian ex- 
pansions for spinning holes in quasi-circular orbit that 
were kindly provided by L. Kidder. The expressions for 
arbitrary spins and masses are lengthy. If one restricts 
attention to equal-mass holes, Mi = M 2 = M, m = 2M, 
fi = M/2, it turns out that only the sum of the spins 
enters the (post) -Newtonian expansions. In terms of 



S!+S 2 

M 2 : 



(13) 



The "++ sequences" with two co-rotating holes. and with L being the unit-vector parallel to L, the 
I (post) 2 -Newtonian expansions become 



f* = --{m£l) 2 'H 1 - ^(m^) 2/3 + ~(t ■ s)(mO) 



384 



—) = {mn)- 2 ln 1 + 2(L ■ sXmtt) 1 / 3 + ( — + s 2 ) (mfi) 2 / 3 
fim J \12 



3(L • s) 2 - s 2 (mfi) 



N4/3 



+ i(L • s)(mfi) |(L • s) 2 - Is 2 ) (mn)^ j. 



(14a) 
(14b) 
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These expressions are plotted in Figs. |l|-|3| together with 
our results from the effective potential method. There is 
remarkable agreement. 

The sum S1 + S2 is zero for all H — sequences with equal 
spin magnitudes, so (post) 2 -Newtonian theory predicts 
that the H — sequences are identical to the non-rotating 
sequence. This is remarkable, and indeed, in Figs. 
the H — sequences are close to the ++0.0 sequence. How- 
ever, a closer look reveals a systematic behavior from 
which we can gain some insight into our assumptions. 
For fixed angular momentum J/^tm, consider the differ- 
ence in binding energy between a point on a H — sequence 
and a point on the non-rotating 0.0 sequence, 

AE b /(i(S) = —(+-S) - ^(0). (15) 

In Fig. |], AEb//j,(S) is plotted as a function of spin 
for several values of angular momentum J/fim. AEb/fi 
varies as the fourth power of spin. This might be a phys- 
ical effect beyond (post) -Newtonian expansions, but 
for the following reason it seems likely that one of our 
assumptions introduces a non-physical contribution to 
AEb/fii too. Figure H strongly suggests that AEb/ (J, is 
converging to a non-zero value as J//xm (and thus sep- 
aration) increase, indicating that there is a contribution 
to AEf, that is independent of the separation of the holes. 
For all spin configurations, Eb must approach zero in the 
limit of large separation, therefore any physical contri- 
bution to AEb should decrease with separation. More- 
over, a coupling between the holes, physical or unphysi- 
cal, will give rise to a separation-dependent contribution 
to AEb/ 11. Therefore the separation-independent contri- 
bution must be a non-physical effect due to properties of 
each isolated hole. A likely candidate is the underlying 
assumption of conformal flatness. At large separations 
each hole should resemble a Kerr black hole, which is not 
conformally flat. 

Since the Kerr metric is the unique stationary state for 
a spinning black hole, if the conformally flat initial data 
for a single hole were evolved, the metric would relax to 
the Kerr metric and emit some gravitational radiation. 
Therefore the total energy contained in our initial data 
slices is larger than in a more faithful conformally non- 
flat data slice and AEb/ n should be positive, which it 
indeed is. 

We conclude that at large separations AEb is contam- 
inated by an unphysical contribution because of the con- 
formal flatness assumption. At small separation there 
might be additional physical contributions beyond the 
(post) 2 -Newtonian order. 

B. Behavior at small separations — ISCO 

In this section we report the key results of this work 
- the spin dependence of the innermost stable circular 
orbit. As we will see, the interpretation of our data at 



small separations is somewhat complicated. At large sep- 
arations, the assumptions and approximations we have 
used are reasonable, except for the assumption of con- 
formal flatness when the holes are spinning. At small 
separations, the interaction between the two black holes 
becomes relatively strong, and our approximations begin 
to break down. Near the ISCO, we must evaluate the 
quality of our assumptions to determine how reliable our 
results are. 

In the neighborhood of each tentative ISCO, we com- 
pute a set of complete effective potential contours. These 
are shown in Fig. ^[ In each plot, the binding energy 
Eb/fi is shown as a function of separation t/m for sev- 
eral different values of angular momentum J/nm. Also 
plotted is the sequence of quasi-circular orbits passing 
through the minima of the effective potential. Figure o 
shows the non-rotating sequence ++0.0, one example 
each of a — and a H — sequence, and three ++ sequences 
with different spin magnitudes. 

Examining the constant J contours of the effective po- 
tential for fixed spin configurations, we find that they 
fall into three regimes separated by critical values that 
we will label J a and Jb- Contours with J > J a ex- 
hibit a single minimum positioned at large separation 
£/m. This minimum moves inward as the angular mo- 
mentum decreases, i.e. the holes approach each other as 
angular momentum and energy are radiated away. We 
call this the "outer" minimum. As J passes through the 
critical value J a, a new "inner" minimum appears in- 
side the outer minimum. In this region, contours of the 
effective potential have two minima separated by a lo- 
cal maximum. The maximum corresponds to the well 
known unstable circular orbit of a Schwarzschild black 
hole. As J decreases further, J a > J > Jb, the maxi- 
mum moves outward whereas the outer minimum contin- 
ues to move inward - the quasi-circular orbit associated 
with the outer minimum continues to shrink. As J passes 
through the second critical value Jb the outer minimum 
and the maximum meet in an inflection point and disap- 
pear. The quasi-circular orbit associated with the outer 
minimum disappears and this inflection point is identi- 
fied with the ISCO. For J < Jb, only the inner minimum 
remains. 

This behavior for the non-rotating sequence was al- 
ready found in [Q. There, the inner minimum was dis- 
missed as unphysical, since the underlying assumptions 
become weaker at small separations of the holes, and 
since a common event horizon might form. We will dis- 
cuss this "unphysical" region and the possibility and con- 
sequences of the formation of a common event horizon 
below. But first we continue discussing the behavior of 
the effective potential for different spin configurations. 

As we increase the spin magnitude for the — config- 
urations, the two critical angular momentum values J a 
and Jb move away from each other. We see a more pro- 
nounced local maximum and the Ef, curves look similar 
to the effective potential of Schwarzschild for a larger in- 
terval of angular momenta. The ISCO moves outward to 
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larger separations as spin increases. 

Conversely, as we increase the spin magnitude for 
the ++ configurations the interval (Jb,Ja), where two 
minima and a local maximum exist becomes smaller. 
Slightly above S/M 2 = 0.17, J a and Jb merge and for 
S/M 2 > 0.17, the regime with two minima and a maxi- 
mum is not present. Figure^ illustrates the small interval 
(Jb,Ja) with an enlargement of the ++0.17 sequence. 
As long as the regime with two minima and a maximum 
is present, we can still define the ISCO by the inflection 
point. It moves towards smaller separation of the holes as 
the spin is increased. However, since the inflection point 
ceases to exist at some spin magnitude, we cannot define 
an ISCO for all S/M 2 . Therefore the ++ sequences dis- 
played in Figs. [l]-|^ do not terminate. Furthermore, we 
need a more careful analysis to determine whether the 
ISCO properties for spin magnitudes close to the critical 
value S/M 2 w 0.17 are reliable. 

The H — configurations are very similar to the non- 
rotating one. Given the weak dependence on spin within 
the H — sequences, this is not surprising. We do not con- 
sider the H — configurations further. 

Figure ^ and Table | summarize the orbital parameters 
at the ISCO as a function of spin for the — sequences 
and the ++ sequences. The numerical errors in Eb/fi, 
L/ ' [im and J / \im are less than 1 per cent, while m£l and 
i/m are accurate to a few percent. However, for the ++ 
sequences the systematic errors of our approach might be 
much larger. The table also includes ISCO parameters 
for a test mass orbiting a Kerr black hole obtained from 
formulas in p0| . 



C. Common apparent horizons 

A common event horizon might be responsible for the 
strange behavior of the effective potential at small sep- 
arations, because once a common event horizon forms, 
there are no longer two distinct black holes. It would be 
helpful to know the critical separation where a common 
event horizon first forms. However, in order to locate the 
event horizon, knowledge of the complete spacetime is 
needed. In the present case, only data on one time-slice 
is available, and so we can only search for common ap- 
parent horizons. Since the event horizon must lie outside 
the apparent horizon, the formation of a common appar- 
ent horizon places a firm bound on the formation of an 
event horizon. 

Searches for a common apparent horizon were carried 
out for several spin configurations. Details of the appar- 
ent horizon finder and the method used to discern the 
formation of a common apparent horizon are given in 
the Appendix. In Table [n], the results of the apparent 
horizon searches are listed. 

For fixed spin configurations the common apparent 
horizon forms at larger separation for larger angular mo- 
mentum. This can be seen from the — 0.25 and +- 1-0.0 



sequences. For varying spins and angular momentum 
close to the ISCO values, the proper separation between 
the throats at the formation of the common apparent 
horizon depends weakly on the spin. It decreases from 
i/m « 2.3 for the — 0.37 sequence down to i/m s» 2.0 
for the ++0.17 sequence. 

Notice that the segment of parameter space where com- 
mon apparent horizons form does not include the se- 
quence of quasi- circular orbit configurations. Indeed, the 
common apparent horizons form at a separation inside 
the inner minimum where the effective potential increases 
with decreasing separation. 

The search for the onset of common apparent horizons 
also provides the actual surfaces. In Fig. || some appar- 
ent horizon surfaces just inside the formation of a com- 
mon apparent horizon are plotted. The circles represent 
the throats of the holes. The solid lines represent a cut 
through the plane of orbital motion of the holes, arrows 
indicating the direction of linear momentum of the holes. 
The dashed lines are cuts through the plane perpendic- 
ular to the plane of motion and parallel to the spins of 
the holes. We find that the apparent horizons lag behind 
the orbital motion, with the amount of lag being larger 
for counter-rotating than for co-rotating holes. 

IV. DISCUSSION 

We found that the effective potential contours at very 
small separation increase with decreasing separation. 
This is in contrast to the usual shape of the effective 
potential for a Schwarzschild or a Kerr black hole, which 
tends to — oo at sufficiently small separations. 

This behavior can be interpreted in the light of the 
common apparent horizon searches. The common ap- 
parent horizon that was found to form at a small sepa- 
ration of the holes might influence the observed effective 
potential as follows: The event horizon must lie outside 
the apparent horizon. Therefore a common event horizon 
must form before a common apparent horizon forms. To 
accomplish this the event horizons around the individ- 
ual holes must grow towards this common event horizon. 
Thus, even before formation of a common event horizon, 
the individual event horizons will no longer be close to the 
individual apparent horizons and the areas of the event 
horizons of the individual holes must be larger than the 
areas of their apparent horizons. Therefore, equations 
and (||) will imder-estimate the mass of the holes. 
We denote this underestimate by AM. Consider the ef- 
fect this underestimate of M has on the binding energy. 
The numerator of (JsJ) will be cwer-estimated by a relative 
amount of 

2AM _ 4 AM AM 
\E ADM - 2M\ ~ \E b /fi\ M > ~M~' " (>! 

At the same time, the denominator of (Q) and the de- 
nominator of the scaled angular momentum ( f7c| ) change 
too, leading to an underestimate of the binding energy 
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EbJ '/J,. However, the relative changes of these denomina- 
tors are only of the order of AM/M, so that the overes- 
timate from Eqn. (|lFJ) dominates. It might well be that 
this overestimate is so large that it counter-balances the 
decreasing effective potential that one might expect in 
analogy to Schwarzschild or Kerr black holes. 

This idea leads to the following picture to explain the 
observed effective potential curves: At large separation 
of the holes, the masses of the holes and the effective po- 
tential are reliable and we see an effective potential that 
looks similar to a Schwarzschild black hole. Consider, 
for example, the ++0.0 sequence: For J slightly above 
its ISCO value we see the (outer) minimum of the stable 
quasi-circular orbit and a maximum corresponding to an 
unstable circular orbit. As J increases, the stable cir- 
cular orbit moves outwards and the unstable one moves 
inwards. Once the maximum corresponding to the un- 
stable orbit moves too far in, the AM/M contamination 
of the effective potential "eats up" the maximum and it 
disappears. 

Now we turn on spin. We found that a common appar- 
ent horizon forms at approximately the same proper sepa- 
ration, independent of the spin of the holes. It seems rea- 
sonable that the AM/M error is also weakly dependent 
on the spin, and also the separation of the holes, where 
AM/M becomes significant. For the — sequences the 
ISCO moves to larger separations. Thus the maximum 
in the effective potential (the unstable orbit) will survive 
for a larger range of separations and angular momenta 
J. Conversely, for the ++ sequences, the ISCO moves 
inwards, closer to the separation where AM/M becomes 
significant. The maximum in E^/fi is lost almost imme- 
diately, and in the extreme limit of S/M 2 > 0.17, it does 
not show up at all. 

This scenario is sufficient to capture the complete be- 
havior of the effective potential as a function of J and 
spin. What does this picture imply for the validity of 
our ISCO results from Table |? We expect that AM/M 
decays rapidly with increasing separation, so the ISCO 
data for the non-rotating sequence ++0.0 as well as for 
the — sequences should be sound. However, because 
AM changes the characteristic behavior for the ++ con- 
figurations even for S/M 2 < 0.17, the ++ sequences will 
be affected. Let us consider how these changes affect our 
estimates of circular orbits. 

Figure ^| illustrates the effect of the AM/M contami- 
nation on the effective potential contours. As we noted 
above, the AM/M contaminations of the binding energy 
overestimates the binding energy of an effective potential 
contour. Since this error increases as the separation de- 
creases, our estimates for the separation at a given value 
of angular momentum are also too high, and our esti- 
mates of the orbital angular velocity mQ are too low. 
Unfortunately, we cannot determine whether our esti- 
mates for the location of the ISCO are too high or too 
low. While our estimates for the separation of a given 
orbit are too high, we see that the true ISCO will occur 
at a larger value of the total angular momentum than we 



estimate. These effects oppose each other. 

The angular momentum at the ISCO, J / fim, increases 
with spin for the ++ configurations. It is interesting to 
examine whether the final black hole resulting from a 
merger of such a spinning binary black hole can violate 
the Kerr limit on spin of a black hole. From (^) we find 



M = 



M 2 



1 - 



S 2 

M 4 



(17) 



By the area theorem, the final irreducible mass must 
satisfy M 2 r j > 2M 2 r , where equal mass holes were as- 
sumed. The final angular momentum cannot exceed the 
angular momentum at the ISCO, J/ < J. With these two 
constraints and by virtue of the Christoudoulou formula 
(0), we find 



M 2 



< i 



(J/fim)" 



4 1 



sjl - (S/M 2 ) 2 



(18) 



A Kerr black hole has always M 2 /M 2 r < 2 with equal- 
ity in the extreme Kerr limit. With data from Table || 
we find for the — 0.50 sequence M 2 /M 2 rf < 1.43 and 
for the ++0.17 sequence M 2 f /M 2 rJ < 1.61. These val- 
ues correspond to spin parameters of J/M 2 < 0.92 and 
J/M 2 < 0.97, respectively. Hence the merged black hole 
might be close to the Kerr limit, but will not violate it. 



A. + Sequences and conformal flatness 

The (spin) 4 effect illustrated in Fig. || suggests that the 
assumption of conformal flatness might lead to inaccu- 
rate results. This is particularly important for analysis 
of gravitational waves. As seen in Fig. |], for spinning 
holes with S/M 2 ~ 0.50 the assumption of conformal 
flatness results in an unphysical gravitational wave con- 
tent of the order of — 2 • 10~ 3 p ~ 5 • 10~ 4 m. This is 
less than 0.1 percent of the total mass and a few percent 
of the binding energy Et>. If the gravitational energy ra- 
diated away is less than 1% of the total mass, then the 
gravitational wave content due to an unsuitable initial 
data slice is a significant contamination. 



V. CONCLUSION 

In this work, we have constructed sequences of quasi- 
circular orbits for equal-sized, spinning black holes. At 
large separations, the results we have obtained match 
well with (post) 2 -Newtonian expansions, although there 
is a clear contamination of the data because of the as- 
sumption of conformal flatness. The main results of this 
paper, displayed in Tabic | and Fig. [?], reveal the be- 
havior of the ISCO for the cases where the spins of the 
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holes are either both co-rotating (++) or counter-rotating 
( — ) with respect to the orbital motion. For co-rotation, 
the ISCO moves inwards with increasing spin and the or- 
bital angular frequency increases. For counter-rotation 
the ISCO moves outward and the orbital angular fre- 
quency decreases. In fact, we find that the orbital an- 
gular frequency changes by almost a factor of 2 between 
the — 0.50 sequence and the ++0.08 sequence. We have 
noted a systematic error in our results that has its origins 
in an underestimation of the mass of each black hole when 
they are close together. For the ISCO, this implies that 
our results are most accurate (ignoring the errors due 
to conformal flatness) when the holes have large counter- 
rotating spins, and the error increases as we move to con- 
figurations with large co-rotating spins. In fact, the error 
becomes so large in the +- 1- sequences that our method 
cannot locate the ISCO when S/M 2 > 0.17. 

Our results clearly show the need to give up the sim- 
plifying assumption of conformal flatness if we are to 
construct astrophysically realistic black hole initial data. 
This is certainly not a new realization, but this is the first 
time that the effects of the conformal flatness assumption 
have been seen so clearly in the context of black hole bi- 
naries. Work toward more astrophysically realistic ini- 
tial data has begun [^l) . This improvement in the initial 
data is needed for all separations. It remains to be seen 
what impact this improvement will have on the process 
of locating quasi-circular orbits when the holes are close 
together. It is likely that the systematic underestimate 
of the mass will still be significant. If so, an improved 
method for locating quasi-circular orbits and the ISCO 
will be useful. 
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APPENDIX: COMMON APPARENT HORIZONS 

Here we provide details of the apparent horizon (AH) 
finder. We use the AH finder described in jlifl . The AH 
surface is expanded in spherical harmonics up to some 
order L. The apparent horizon, as a marginally outer 
trapped surface, has everywhere vanishing expansion and 
is located by minimizing the square of the expansion over 
the surface. We use convergence with increasing expan- 
sion order L to diagnose the formation of a common AH. 
Therefore high expansion orders L are needed as well as 



reliable convergence of the minimization routine to the 
true minimum of the square of the expansion. 

The Powell minimization used in ]l2] is too slow for 
high-order expansions. We replaced it by a DFP method 
with finite difference approximations of the Jacobian H . 
For the modest expansion order L — 6, DFP is already 
ten times faster than Powell's method. 

Furthermore, we take advantage of the symmetries of 
the AH surface. The holes are located along the z-axis 
at z = ±/3/2. Their linear momenta point in the in- 
direction and the spins are directed along the ±y-axis. 
It is straightforward to show that these choices imply 
that the AH surface is invariant under reflection at the 
xz-plane, y — » —y. This symmetry constrains the coeffi- 
cients Aim of the expansion in spherical harmonics to be 
real. Moreover, for the ++ and — configurations with 
equal sized holes and equal spin magnitudes, the config- 
uration is symmetric under rotation by 180° around the 
y-axis, this is (x,y,z) — > (—x,y,—z). Both symmetries 
together force Ai m — for odd I and Ai m to be real for 
even I. Hence the number of free parameters in the min- 
imization routine can be reduced by almost a factor of 
four. 

To prevent convergence to spurious local minima, it is 
vital that the function that is minimized be as smooth 
as possible. Therefore we use second order spline inter- 
polation to provide the required data for the AH finder. 
Compared to bicubic interpolation, the spline interpola- 
tion somewhat decreased the number of iterations needed 
in the minimization routine, but more importantly it sig- 
nificantly reduced the probability of getting stuck in a 
local minimum. In addition, many rays were used to 
reduce the anisotropies introduced by the discrete po- 
sition of the rays. Finally, we distribute the rays non- 
uniformly in solid angle. The reason for this is simple: 
The common AH surface will be very oblate along the 
z-axis, since it must encompass the two throats located 
along the z-axis. The polar regions of the AH surface 
are close to the throats and the conformal factor changes 
rapidly. These regions are particular important, but the 
standard distribution uniform in cos 9 places relatively 
few rays in the polar regions. Therefore we implemented 
a procedure that distributes the rays in proportion to an 
arbitrary ray-density function f(0). A uniform distribu- 
tion of rays is represented by f(9) = const., whereas we 
used f(9) = 1 + cos 2 9, resulting in a doubled density of 
rays close to the poles. 

With the improved AH finder, we performed extensive 
tests with various numbers of rays. As a rule of thumb, 
about ten times more rays as free minimization param- 
eters are necessary to ensure reliable convergence to the 
true minimum of the square-sum of the expansion. 

We used expansions up to order L = 16 and up to 
64x48 rays (64 in 8 direction, 48 in 0). We perform a 
set of AH searches, starting at L = 2 and increasing L 
by 2 between searches. The result of the previous search 
is used as the initial guess for the next higher expansion 
order. Such a set of expansions from L = 2 to L = 16 
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takes typically about 2 hours CPU time on a RS6000 
processor. 

A disadvantage of an AH finder based on a minimiza- 
tion routine is that the minimization routine will always 
find a minimum. It does not matter whether there ac- 
tually is a "true" apparent horizon, or whether there is 
only a surface with small but non-zero expansion. And 
even for a true AH, the result of the minimization will 
be non-zero because of the finite grid resolution in the 
underlying elliptic solver and finite expansion order in 
spherical harmonics. Therefore we need a method to dis- 
cern a "true" AH from a mere minimum in the square of 
the expansion. 

For a true AH, the square of the expansion is exactly 
zero, therefore we expect that the residual of the min- 
imization tends to zero as the resolution of the elliptic 
solver and the expansion order L are increased. With 
increasing L, the error in the approximation of the sur- 
face by spherical harmonics should decrease exponen- 
tially. On the other hand, for a mere minimum in the 
expansion, we expect that the residual of the minimiza- 
tion tends towards a non-zero limit as the resolution of 
the elliptic solver and the expansion order L is increased. 
We use this signature to discern the formation of a com- 
mon apparent horizon. 

Figure ^ shows the residual of the minimization for 
various values of L and different separations (3. The solid 
lines represent configurations at different separations of 
the holes. They are labeled by the background separation 
of the holes, (3. Each solid line represents a set of min- 
imizations with varying expansion order 2 < L < 16 on 
the same initial data set. At large separations, (3 > 4.5, 
the residual of the minimization becomes independent of 
L for large L. At small separation, (3 — 4.4, the residual 
decreases exponentially through all computed expansion 
orders up to L = 16 - a common AH has formed. 

Neither reducing the number of rays, nor decreasing 



the resolution of the Hamiltonian solver changes the con- 
vergence behavior significantly. This is illustrated by 
some examples in Fig. We conclude that for this 
particular example a common AH first forms between 
(3 = 4.4 and (3 = 4.5. 

Expansions to high order in L are essential for discern- 
ing the formation of a common AH. If one had Fig. [l(] 
only up to expansions up to L — 8, it would be impos- 
sible to decide where the common AH first forms. One 
would probably conclude that the common AH forms at 
larger separations than it actually does. 
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TABLE I. Orbital parameters of the innermost stable circular orbit for equal-mass spinning 
holes. The second through sixth columns give the data obtained in this work, the three columns to 
the right give the data for a test mass orbiting a Kerr black hole. The results for the +- 1-0.08 and 
+- 1-0.17 sequences will have larger systematic errors than the other cases (see text). 



sequence 


t/m 




mO, 


J//j,m 


Lj lira 


E b /fi 


L/ firn 


mO 


0.50 


7.05 


-0.0628 


0.100 


2.438 


3.438 


-0.04514 


3.8842 


0.04935 


0.37 


6.68 


-0.0687 


0.107 


2.595 


3.335 


-0.04767 


3.7834 


0.05319 


0.25 


6.17 


-0.0743 


0.120 


2.730 


3.230 


-0.05032 


3.6856 


0.05727 


0.12 


5.58 


-0.0815 


0.139 


2.865 


3.105 


-0.05363 


3.5738 


0.06242 


++0.0 


4.94 


-0.0901 


0.166 


2.976 


2.976 


-0.05719 


3.4641 


0.06804 


++0.08 


4.59 


-0.0975 


0.186 


3.042 


2.882 


-0.05991 


3.3870 


0.07237 


++0.17 


3.93 


-0.1087 


0.235 


3.103 


2.763 


-0.06337 


3.2957 


0.07793 



TABLE II. Summary of the common apparent horizon 
searches. Listed are the sequences and values of orbital an- 
gular momentum for which an apparent horizon search was 
carried out. The apparent horizon was found to form at a 
separation l\/m < £/m < lijm. 



Sequence 


L/fim 


h/m 


£ 2 /m 


0.37 


3.38 


2.32 


2.38 


0.25 


3.10 


2.20 


2.25 


0.25 


3.34 


2.24 


2.29 


++0.0 a 


0.0 


1.89 




++0.0 


2.94 


2.08 


2.13 


++0.0 


3.00 


2.08 


2.13 


++0.08 


2.84 


2.03 


2.08 


++0.08 


2.92 


2.03 


2.08 


++0.17 


2.79 


1.98 


2.03 


++0.25 


2.70 


1.96 


2.01 


+-0.25 


3.00 







a From ]13j, which found a critical separation (3 — 4.17. This 
corresponds to a proper separation of tjm m 1.89. 
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FIG. 1. Sequences of quasi-circular orbits for different 
spin configurations. Plotted is the binding energy Eb/n vs. 
the angular momentum J '/ ' jxm along the sequences. The solid 
lines represent the data, the dashed lines are the results based 
on (post) 2 -Newtonian theory. As discussed later in this paper, 
the effective potential method could not locate an ISCO for 
the ++0.25 and ++0.50 sequences, although we believe each 
sequence should terminate in one. 
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FIG. 3. Sequences of quasi-circular orbits for different 
spin configurations. Plotted is the angular momentum J/fim 
vs. the orbital angular frequency mfi along the sequences. 
The solid lines represent the data, the dashed lines are the 
results based on (post) 2 -Newtonian theory. As discussed later 
in this paper, the effective potential method could not locate 
an ISCO for the ++0.25 and ++0.50 sequences, although we 
believe each sequence should terminate in one. 
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FIG. 2. Sequences of quasi-circular orbits for different spin 
configurations. Plotted is the binding energy E b //j, vs. the 
orbital angular frequency tnfi along the sequences. The solid 
lines represent the data, the dashed lines are the results based 
on (post) 2 -Newtonian theory. As discussed later in this paper, 
the effective potential method could not locate an ISCO for 
the ++0.25 and ++0.50 sequences, although we believe each 
sequence should terminate in one. 
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FIG. 4. Difference in binding energy AE b /fi between H — 
sequences and non-rotating sequence as a function of spin of 
the +— sequence for fixed angular momentum J//xm. Each 
curve is labeled by its value of J/fim. J/fim = 3.01 is very 
close to the ISCOs that have J/^im ~ 2.98. J/pim = 3.65, 
3.35, 3.15 and 3.01 correspond to a separation of £/m w 12.3, 
9.6, 7.7 and 6.1, respectively. 
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S 4 6 8 2-168 

t/m l / m 

FIG. 5. Constant J / \im contours of the effective potential 
Eb/fi as a function of separation t/m for various spin config- 
urations. The curves are spaced in steps of AJ/pm = 0.02 

except for the 0.25 and the H — h0.17 configurations, which 

have steps of 0.04 and 0.01, respectively. Also plotted is the 
sequence of quasi-circular orbits connecting the minima of the 
effective potential. 



-0.5 -0.25 

S y ■ L — S p ' L 

FIG. 7. Values of several physical parameters at the ISCO 
of the ++ and sequences. Plotted are the binding en- 
ergy Et/fi, the orbital angular frequency mQ, the total angu- 
lar momentum J/fj,m and the proper separation between the 
holes, t/m as a function of spin S/M 2 on the holes. The ++ 
sequences are plotted along the positive part of the horizontal 

axis, the sequences along the negative part as —S/M 2 . 

The vertical axes on the left side belong to E^/jj, and J / pan. 
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FIG. 6. Enlargement of the ++0.17 sequence of Fig. || 
The displayed effective potential contours (top to bottom) 
correspond to angular momenta J / prni — 3.12, 3.11, 3.104, 
3.103, 3.102, 3.10, 3.09 and 3.08. Also shown is the sequence 
of quasi-circular orbits. 
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FIG. 8. Shapes of the common apparent horizons for dif- 
ferent spin configurations. Circles denote the throats of the 
holes. The solid lines are cuts in the plane of orbital motion 
(arrows indicating the direction of motion), the dashed lines 
represent cuts normal to the plane of motion. 
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r/M 

FIG. 9. Illustration of the effects of a systematic under- 
estimation of Eb/fi. The dashed lines represent the observed 
effective potential contours for some values of J. The points 
A, B, and C correspond to circular orbits. The ISCO is at D. 
Assuming that the true binding energy is smaller, with the 
deviation increasing as the separation decreases, yields true 
effective potential contours similar to the solid lines. The 
true circular orbits are at E and F and the true ISCO is at 
G. We find that the minima of the true contours will lie at 
smaller separation (for the same J). The angular frequency 
is given by Q — dE b /dJ. Using the points A and B, we see 
that the observed dEb is smaller than the true one, so we un- 
derestimate mfl. For fixed J, true circular orbits will occur 
at smaller separation, but the true ISCO will appear at larger 
J than we have observed. These effects counteract each other, 
making it impossible to predict their effect on the true ISCO. 
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FIG. 10. Residual of the minimization in the AH finder as 
a function of expansion order L. The number of rays used was 
Ne — 64, Ntf, — 48. The different solid lines represent different 
separations of the holes along an effective potential contour 
with J I jj,m — 3.29 on the H — h0.25 sequence. The dashed lines 
are the results of minimizations with Ng = 48, N$ = 32. The 
dotted-dashed lines show examples of minimizations at lower 
grid resolution and Ne =64, Nj, = 48. 
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